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Magnetowave induced plasma wakefield accel- 
eration (MPWA) in a relativistic astrophysical 
outflow has been proposed as a viable mechanism 
for the acceleration of cosmic particles to ultra 
high energies [1]. Here we present simulation 
results that clearly demonstrate the viability of 
this mechanism for the first time. We invoke the 
high frequency and high speed whistler mode 
for the driving pulse. The plasma wakefield 
so induced validates precisely the theoretical 
prediction. We show that under appropriate 
conditions, the plasma wakefield maintains very 
high coherence and can sustain high-gradient ac- 
celeration over a macroscopic distance. Invoking 
gamma ray burst (GRB) as the source, we show 
that MPWA production of ultra high energy 
cosmic rays (UHECR) beyond ZeV (10 21 eV) is 
possible. 

The origin of ultra high energy cosmic rays (UHECR) 
is a long-standing mystery in astrophysics. Thus far, the 
theories that attempt to explain the origin of UHECR 
can be broadly categorized into the "top-down" and 
the "bottom-up" scenarios. Each scenario faces its own 
theoretical and observational challenges 2]. Precision 
measurements [H, H| on the yield of air-shower induced 
fluorescence lend support to the energy calibration of 
the HiRes observations Recent data from the Auger 
Observatory Q, also fluorescence normalized, exhibit a 
similar location of the "ankle" and super-GZK steepening 
as HiRes. The lack of a strong super-GZK flux reduces 
the need for top-down exotic models. If these models 
are indeed disfavored, then the challenge to find a viable 
"bottom up" mechanism to accelerate ordinary particles 
to and beyond 10 20 eV becomes more acute. 

Shocks, unipolar inductors and magnetic flares are the 
three most potent, observed, "conventional" accelerators 
that can be extended to account for ~ ZeV(= 10 21 eV) 
energy cosmic raysQ- Radio jet termination shocks and 
gamma ray bursts (GRB) have been invoked as sites 
for the shock acceleration, while dormant galactic cen- 



ter black holes and magnetars have been proposed as 
sites for the unipolar inductor acceleration and the flare 
acceleration, respectively. Each of these models, how- 
ever, presents problems 7]. Evidently, novel acceleration 
mechanisms that can avoid the difficulties faced by these 
conventional models should not be overlooked. 

Plasma wakefield accelerators H Q are known to pos- 
sess two salient features: (1) The energy gain per unit 
distance does not depend (inversely) on the particle's in- 
stantaneous energy or momentum. This is essential to 
avoid the gradual decrease of efficiency in reaching ul- 
tra high energies; (2) The acceleration is linear. Bend- 
ing of the trajectory is not a necessary condition for 
this mechanism. This helps to minimize inherent energy 
loss which would be severe at ultra high energy. How- 
ever, high-intensity, ultra-short photon or particle beam 
pulses that excite the laboratory plasma wakefields are 
not readily available in the astrophysical setting. It was, 
however, proposed[lj that large amplitude plasma wake- 
fields can instead be excited by the astrophysically more 
abundant plasma "magnetowaves" , whose field compo- 
nents are magnetic in nature (\B\ > \E\). Protons can be 
accelerated beyond ZeV energy by riding on such wake- 
fields. Though attractive, this concept has never been 
validated through self-consistent computer simulations. 
In this Letter, we report on the plasma particle-in-cell 
(PIC) simulation results that confirm the magnetowave- 
induced plasma wakefield acceleration (MPWA) concept 
for the first time. 

Magnetized plasmas support a variety of wave modes 
propagating at arbitrary angles to the imposed magnetic 
field. For the specific case of wave modes propagating 
parallel to the external magnetic field, the electromag- 
netic waves become circularly polarized and the disper- 
sion relation is [Tol 
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where the upper (lower) signs denote the right-hand (left- 
hand) circularly polarized waves. tu p — y/ '4ire 2 n p /m e 
is the electron plasma frequency, uj c — eB/m e c is the 
electron cyclotron frequency and the subscript i denotes 
the ion species. Each polarization has two real solutions 
with high and low frequency branches and both have a 
frequency cutoff which forms a forbidden gap for wave 
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propagation. The right-hand polarized, low frequency 
solution is called the whistler wave which propagates at 
a phase velocity less than the speed of light. When the 
magnetic field is sufficiently strong such that lu c 3> u> p , 
the dispersion of the whistler mode becomes more linear 
over a wider range of wavenumbers with phase velocity 
approaching the speed of light (see Fig.l). The E and B 
components of the wave are now comparable in strength. 
In this regime the travelling wave pulses can maintain 
their shape over macroscopic distance, a condition desir- 
able for plasma wakefield acceleration. 

The ponderomotive force in a magnetized plasma has 
been well studied [ll|. Applying the dispersion relation 
for the whistler wave, with the ion motion neglected, we 
obtain the ponderomotive force acting on an individual 
electron as 
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where E w (£) is the amplitude of the whistler wave- 
packet, and Q = z — v g t the co-moving coordinate for 
the driving pulse. Note that E w is perpendicular to z. 

Combining this equation with the continuity equation 
and the Poisson equation, the longitudinal electric field 
in the plasma, i.e., the plasma wakefield, can be solved 
and it reads 
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where E w is the maximum value of E n ,((). An expres- 
sion similar to Eq. ([3|) has been obtained for laser-induced 
wakefield in a magnetized plasma fl2l|. For a Gaussian 
driving pulse with E W {C) = E w exp(— ( 2 /2<r 2 ), it can be 
shown that behind the driving pulse, i. e., |£| 3> er, 
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It is customary to express the plasma wakefield in terms 
of the Lorentz invariant "strength parameter" of the driv- 
ing pulse, ao = eE w /m e ctu, and the "wavebreaking" 
field, E w b = m e cu!p/e. Assuming the driving pulse fre- 
quency is centered around u> and its speed v g w oj/k, the 
maximum wakefield, or the acceleration gradient, attain- 
able behind the driving pulse is then 
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Relative to the conventional wakefields, G is enhanced 
by a factor k 2 c 2 /(u> — uj c ) 2 when w approaches u c . The 
above expression is derived under the assumption of lin- 
ear plasma perturbation, i.e., ao < 1. Research made 



FIG. 1: (a) Frequency and (b) phase velocity versus wavenum- 
ber for different magnetic field strengths. The vertical solid 
line is the mean value of the pulse wavenumber that was cho- 
sen for the PIC simulation, and the dashed lines its range. 



over the past two decades in plasma wakefields has firmly 
established the generalized wakefield amplitude for all 
values of ao [H[ , and adapting its form we obtain, 
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Note that in the nonlinear regime (ao 3> 1) the wakefield 
is no longer sinusoidal in £ but saw-tooth-like. 

We have conducted computer simulations to study the 
MPWA process driven by a Gaussian driving whistler 
pulse described above. Our simulation model integrates 
the relativistic Newton-Lorentz equations of motion in 
the self-consistent electric and magnetic fields determined 
by the solution to Maxwell's equations [3, [3- The 
4-dimensional phase space (z,p x ,p y ,p z ) is used for the 
charged particle dynamics and a uniform external mag- 
netic field, -Bo, is imposed in the z-direction. In order 
to sustain the driving pulse shape propagating in the 
plasma, it is important to have modes of the pulse trav- 
elling with similar phase velocities. Therefore, we used 
a wavepacket with Gaussian width a = 80 A/ v^, where 
A is the cell size taken to be unity, and the wavenum- 
ber k = 27r/60A. The normalized physical parame- 
ters ll) c /u> p = 6 and mi/m e — 2000 were taken and 
for a uniform background plasma with electron collision- 
less skin depth, c/lu p A — 30, this gives uj/lj p = 2.98 
and Vg/c ~ w/ck = 0.95. Other numerical parame- 
ters used are: total number of cells in the z-direction, 
L z = 8192A = 273c/o7 p , average number of particles per 
cell was 10, and the time step oj p At = 0.1. The fields 
were normalized by (l/30)E w j ) . 

We set the maximum amplitude E w — 10, which 
gives the normalized vector potential ao = eE w /m e cu> = 
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FIG. 2: A snapshot of the plasma wakefield induced by the 
whistler pulse. E x is in blue and E z in red. 



FIG. 4: Validation of the functional dependence of G in 
Eq.©, where E w b is in units of (l/30)m e cujp/e. 




FIG. 3: The dispersion relation of the driving pulse from PIC 
simulation. Liner dispersion curves for the L-wave (orange) 
and the whistler wave (pink) are superimposed. 



0.11 <C 1. Thus the wakefield in our simulation is in 
the linear regime. The pulse was initialized at zq — 
500A = 16.66c/cj p . To avoid spurious effects, we grad- 
ually ramped up the driving pulse amplitude until t — 
lOOWp 1 , during which the plasma feedback to the driv- 
ing pulse was ignored. After this time, the driving pulse- 
plasma interaction was tracked self-consistently. As the 
dispersion relation in this regime is not perfectly linear, 
there was a gradual spread of the pulse width. Thus x 
and E w of the driving pulse decrease accordingly. As a 
result, the maximum wakefield amplitude, E z , declined 
in time. Even so, it agrees very well with the theoretical 
maximum of E z ~ 0.266(l/30)i? to b. Fig. 2 is a snapshot 
of E x and E z at t = 230w~ 1 . We note that while the 
driving pulse continues to disperse, the wakefield remains 
extremely coherent. 

We sampled the E x (k) of the pulse, after its initializa- 
tion, at every time step and analyzed it in the frequency 
space. Fig|3]shows the co — k intensity generated from the 



PIC simulation driving pulse power spectrum. It is super- 
imposed with the theoretical curves for the left-handed 
circularly polarized electromagnetic wave (L-wave) and 
the whistler wave dispersion relations. We confirm that 
our driving pulse is indeed a whistler wave. 

Next we validate the functional dependence of the ac- 
celeration gradient given in Eq.©. Fig 0] plots the ac- 
celeration gradient versus xk 2 c 2 E w b/(uu — lo c ) 2 . This is 
performed by varying lo c . Wavepackets with the same 
wavenumber (2-7r/60A) and the same maximum ampli- 
tude (E w — 10) are initialized with different B such 
that uj c /u)p = 6,7,8,9,10. The linear fitting slope is 
0.012, which agrees very well with the expected value of 
Oq. We comment that when oj approaches u> c , i. e., when 
the two are in resonance, the wakefield amplitude can be 
dramatically increased by several orders of magnitude. 
In view of the successful validation by our computer sim- 
ulation of the MPWA in the linear regime, we are safe to 
extrapolate it to the nonlinear regime according to the 
theoretical formula in Eq.©. 

We now relate this mechanism to the issue of UHECR 
acceleration. An earlier attempt has been made[l[ where 
gamma ray bursts (GRB) were invoked as the candi- 
date site for plasma wakefield production of UHECR. 
We follow the same approach here. GRBs are gener- 
ally classified into long (t; ~ 10 — lOOsec) and short 
(t s ~ lsec) bursts. Within seconds (for short bursts), 
about £grb ~ 10 50 erg of energy is released through 
gamma rays. We invoke the Neutron Star-Neutron Star 
coalescence as our working model for GRB. Neutron stars 
are known to be compact (i?NS ~ O(10)km) and carry 
intense surface magnetic fields (-Bns ~ 10 13 G). It is 
generally believed that when neutron stars collide, the 
tremendous release of energy creates a highly relativis- 
tic out-bursting fireball (jets) [16], most likely in the form 
of a plasma. We assume such a jet has an open-angle 
of #grb ~ 0.1 and the initial plasma density in the jet 
^GRB ~ 10 26 cm~ 3 . We further assume that such vio- 
lent collision of intense magnetic fields would create se- 
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quence of strong magnetoshocks, where whistler waves 
are imbedded. In the aftermath of such tremendous im- 
pact, the magnetic field-lines would be temporarily shat- 
tered and reoriented with strong poloidal field lines par- 
allel to the axis of the jet. 

From the previous discussion we see that the MPWA 
is most effective when the driving pulse frequency falls 
between the plasma frequency, u> p , and the electron- 
cyclotron frequency, lo c . Let us verify whether and 
where this condition can be satisfied along the GRB jet. 
First we note that due to the conservation of the mag- 
netic flux, the poloidal magnetic field strength decreases 
as 1/r 2 away from the epicenter of GRB. This means 
lo c = eB/m e c oc (Rns/t) 2 - On the other hand, the con- 
tinuity condition requires that the plasma density de- 
creases as 1/r 2 as well. Therefore, ui p oc i?Ns/ r and the 
cross-over between these two parameters does exist at a 
distance R ~ 100i?NS ~ 1000km. This is the "sweet 
spot" where MPWA is the most effective. 

To estimate the plasma wakefield acceleration gradi- 
ent at this "sweet spot", we first note the EM energy 
density of GRB is E^ rb /4tt = £grb/ct s tt(8 G b,bR) 2 ~ 
10 27 erg/cm 3 . We assume that a fraction, 7] a , of this out- 
burst energy goes into the magnetoshocks. We further 
assume that a fraction, rjt, of the magnetoshocks energy 
lies in the whistler mode. This means E^ ~ iJu^-Egrb- 
Based on this, the associated strength parameter ao is: 
a o = \/VaJ]b(eEGRB) / (rriecuj). At the sweet spot where 
u! c ~ lu ~ fee ~ uip, the factor k 2 c 2 /(lu 2 — lu 2 ) is of the 
order unity. Furthermore, the extremely sharp magne- 
toshock fronts would render the form factor \ a l so of the 
order unity. Assume that ao > 1. Then the acceleration 
gradient boils down to, cf. Eq.©, 

G ~ a a E wb — Jf] a rib ( eE ° RB \ E wb . (8) 
V m e CLU J 

To appreciate what this translates into physical require- 
ments, let us assume that the range of the "sweet spot" 



is 5R ~ 0.1R ~ 10i?NS ~ 100km around R where the 
factor k 2 c 2 /(lu 2 — lu 2 ) is of the order unity. Then in or- 
der for MPWA to be responsible for the production of 
UHECR beyond ZeV (10 21 eV), it is necessary that G ~ 
10 14 eV/cm. In turn, the fractions of GRB energy im- 
parted into the whistler mode have to be rj a ~ rjb ~ 10~ 2 . 

As shown in Ref.[l[, the stochastic encounters of the 
test particle with the random acceleration-deceleration 
phases would result in a inverse-square-law spectrum, 
/(£) oc 1/E 2 . The various additional energy loss mecha- 
nisms, such as few-body collision and synchrotron radia- 
tion, would degrade the power-law index to l/£ 2+f3 , with 
< p < 1. 

Our PIC simulations have confirmed the concept of 
plasma wakefield excited by a magnetowave in a magne- 
tized plasma. Different from the laser and particle beam, 
the magnetowaves are medium waves which cannot exist 
without the plasma. MPWA should thus be of interest 
as a fundamental phenomena in plasma physics and an 
alternative approach to plasma wakefield acceleration. 

As a first step, wc investigated MPWA in the parallel- 
field configuration. Since both poloidal and toroidal field 
components are inevitable in astro-jets, we will further 
investigate plasma wakefield excitation and acceleration 
under the cross-field configuration. In order for MPWA 
to be responsible for the ZeV UHECR production, the 
energy transfer efficiency for GRB is constrained. It 
would be very interesting both observationally and 
theoretically to test whether this constraint is valid. 
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